POP Implementation of Mesoscale eddy diffusivity based on mixing theory

Here are my notes on the development of a parameterization for eddy diffusivity. The scheme is based on the paper by Michael Bates, Ross Tulloch, John Marshall, and Raffaele Ferrari. I'm working with Gokhan Danabasoglu and Matt Long here at NCAR and John Marshall (co-author) at MIT.

The general approach of this effort will be to follow John and Gokhan's implementation notes, isolating each of the terms contributing to the eddy diffusivity given by equation 6 in Bates et al. and comparing the CESM implementation of those terms with Bates, reanalysis, or other observations.

The plots on this web page are averaged from the last 20 years of a 60 year run using the steering parameterization. The experiment was setup similar to the giaf.c112.000 control case and initialized from year 249 of that run. In the standard POP diagnostics below corresponding years were compared between the control and our experiment.

Development Status

As with all new parameterizations getting the code to run and give reasonable results is the first hurdle. The current implementation runs smoothly and the code has been reviewed by Gokhan. Although the ocean state is not in equilibrium yet, the successes and deficiencies of our implementation are evident in this short run.

Good to Go:

Code running and producing reasonable values for many of the diagnostic terms.
CESM Eddy Length scales
First Baroclinic Wavespeed
Zonal mean velocity looks ok compared to ECCO annual average
$(U-c)$ term is close if the phase speed $(c)$ is limited to 20 cm/s

Needs Work:

Zonal Eddy Phase Speed $(c)$
- Excessive from 10S:10N ( 20 cm/s limit imposed )
- Missing westward velocities in ACC which appear in Hughes Data
Eady inverse timescale ($\sigma$)
- Perhaps this is the best we can do given the coarse spatial resolution and daily forcing but this field ultimately feeds in to the derivation of u_rms which is not represented correctly in the region of the western boundary currents.
- The field is too strong and narrow along coastal areas producing very strong eddy diffusivities. We are looking into whether this may be causing the biases in temperature and density gradients in these areas.
u_rms
- As mentioned above this field diverges from obs in the East to West tracks of the Gulf Stream and Kurishio. Due to the way the parameterization is constructed, when u_rms is missing the derived eddy diffusivities are being applied in areas where the zonal velocity is minimum. For the Gulf stream that means a diffusivity being applied at the southern boundary of the Gulf stream.
Suppression Factor
- Spatial patterns seem to match obs for many regions but field too low overall (higher suppression - lower diffusivity values). Supression factor missing in Central Pacific and Atlantic (~40N) at Eastward turning points of WBCs.
Eddy diffusivity
- Spatial patterns seem to match obs for many regions but there are significant errors in Central Pacific and Atlantic (~40N). We are also missing the maximums at 25S and 35S that appear in Bates et al.
Our implementation of the steering parameterization is aggravating biases in the Western Boundary Currents. The resulting Gulf Stream track becomes very zonal and too narrow.

POP standard diagnostics for the new steering case

Time averaged diagnostics for the last 20 years (291_310)

POP timeseries diagnostics for 60 years

New steering case vs giaf control

Implementation Notes on Diffusivity with Steering Level Suppression (Gokhan and John Marshall)

12 November 2014

Starting with

$\quad\quad\quad K=u_{rms}∗L_{mix}$

the general form for diffusivity, $K$, is given by equation (6) of Bates et al. (2014). Namely,

$\quad\quad\quad L_{mix} = {\Gamma * L_{eddy} \over (1 + b1 * |u_{mean} - c|^2 /u_{rms}^2 (z=0)} $

Here,

$u_{rms}$ is the root-mean-square (rms) eddy velocity;

$L_{mix}$ is the mixing length;

$L_{eddy}$ is the eddy diameter (depth independent);

$u_{mean}$ is the mean zonal velocity (resolved);

$c$ is the zonal eddy phase speed (depth independent);

$\Gamma = 0.35$

$b1 = 4$

$\quad\quad\quad\quad\quad\quad$ Eddy Length Scales:

Rossby deformation radius = $L_r = {c_r \over |f|}$

Equatorial Rossby deformation radius $= L_{req} = {\sqrt{c_r \over 2\beta}}$

Rhines scale $= L_{Rh} = {\sqrt{u_{rms} \over \beta}} \sim {\sigma_{vi} \over \beta}$

$c_r$ is the first baroclinic wave speed computed following equation (2.2) of Chelton et al. (1998) with $m=1$;

$f$ is the Coriolis parameter;

$\beta$ is the latitudinal variation of the Coriolis parameter; and

$\sigma_{vi}$ is the Eady growth rate given by

$\sigma_{vi} = {f \over \sqrt{R_i}}$

with $R_i$ the vertically integrated (over the 100 – 2000 m depth range) Richardson number.

So, any one these length scales could be used as an eddy length scale. An alternative is

$L_{eddy} = min (L_r, L_{req}, L_{Rh}).$

$\quad\quad\quad\quad\quad\quad$ Eddy Velocity:

${u_{rms} = alpha*\sigma_{vi}*L_r}$

where $\sigma_{vi}$ is the Eady growth rate based on local Richardson number and $\alpha$ is a scaling constant.

$\quad\quad\quad\quad\quad\quad$ Zonal phase speed:

$c = - \beta * L_r^2$

References:

Bates, M., R. Tulloch, J. Marshall, and R. Ferrari, 2014: Rationalizing the spatial distribution of mesoscale eddy diffusivity in terms of mixing length theory. J. Phys. Oceanogr., 44, 1523-1540, doi: 10.1175/JPO-D-13-0130.1.

Chelton, D. B., R. A. deSzoeke, M. G. Schlax, K. E. Naggar, and N. Siwertz, 1998: Geographical variability of the first baroclinic Rossby radius of deformation. J. Phys. Oceanogr., 28, 433-460.

Tullock, R., J. Marshall, and K. S. Smith, 2009: Interpretation of the propagation of surface altimetric observations in terms of planetary waves and geostropic turbulence. J. Geophys. Res., 114, C02005, doi: 10.1029/2008JC005055.

Questions:

In the 2D implementation, $u_{rms}$ and $u_{mean}$ specifications: upper-ocean vertically or integrated or at $z = 0$?

A: $U_{rms}$ is not depth dependent; both $u_{rms}$ and $u_{mean}$ are for surface only

In the 2D implementation, vertical profile will be specified by $N2(z)$?

A: Yes, $N^2(z) \over N_{ref}(z)$ to be more precise

Local $R_i$ use imbedded in sigma in $u_{rms}$ calculation?

A: No, $u_{rms}$ is depth independent

alpha = ?

A: trial and error

Cancellation of $f$’s in $u_{rms}$ calculation?

Zonal phase speed equation correct? Both $\beta$ and $L_r$ will be positive, producing $c < 0$ always. This appears to be in contrast with Tullock et al. (2009).

A: What is plotted in Bates is (U-c)

$\color{green}{\quad CESM \quad Implementation \quad of \quad Steering \quad Parameterization}$

During the development of the parameterization several tuning modifications were made to the terms listed above to better match the observations.

New implementation of Eady Growth rate (sigma) that avoids undefined values at the equator. See alternate derivation at end of this page.

$\quad\quad\quad \sigma_{vi} = {f \over \sqrt{R_i}}$ is replace by $\sigma_{vi} = \quad {{\require\cancel f m^2} \over \cancel f N}$

Don't use Rhine Scale in calculation of $L_{eddy}$

$\quad\quad\quad L_{eddy} = min (L_r, L_{req}, \cancel{L_{Rh}})$

Eddy phase speed calculation modified to avoid undefined value at the equator and limited to better match Bates paper.

$\quad\quad\quad \cancel{c = - \beta * {L_r^2}}$ is replace with $c = max(- \beta * L_{eddy}^2,-20)$

$u_{rms}$ is hard to capture in this implememtation due to the coarse model resolution and use of mean vs resolved states. We imposed a $5 cm/sec$ minimum on $U_{rms}$ to tune to surface obs of Bates and included a scaling constant $\alpha$ to handle the bias introduced by using mean state instead of resolved state.

$\quad\quad\quad u_{rms} = max(5.,alpha*\sigma_{vi}*L_r)$

$\quad\quad\quad \alpha=4$ The scaling constant for $U_{rms}$

The Eddy diffusivity calculated by our implementation was weak so be bumped up $\Gamma$ from .35 to 1.75 (5x increase)

$\quad\quad\quad \Gamma=.35$ replated with $\Gamma = 1.75$

$\color{green} {Parameterizing \quad the\quad Eddy\quad Length\quad Scale} $

NOTES:

$K=u_{rms}∗{\Gamma * \color{red}{L_{eddy}} \over (1 + b1 * |u_{mean} - c|^2 /u_{rms}^2 (z=0)}$
$L_{eddy} = min (L_r, L_{req}, \require{cancel}\cancel{L_{Rh}})$
Since the two length scales we are using in the parameterization of $L_{eddy}$ depend on the Baroclinic wave speed $c_r$ we will first check $c_r$ against Chelton 1998.

$\quad\quad\quad\quad\quad\quad\quad$ First Baroclinic Wave Speed $c_r$ CESM

In [35]:

Expand Code

$\quad\quad\quad\quad\quad\quad$ Chelton Baroclinic gravity wave phase speed

In [36]:

Expand Code

$\quad\quad\quad$ The barclinic wave speed looks reasonable so lets compare the Paramaterized Rossby Radius to Chelton

$\quad\quad\quad\quad\quad\quad$ First Baroclinic Rossby Radius CESM $\quad (L_{eddy})$

In [37]:

Expand Code

# <!-- collapse=True --> 

fcor = nc.variables['FCOR'][0,:,:]
fcor_nearest = pyresample.kd_tree.resample_nearest(orig_def, fcor, \
        targ_def25, radius_of_influence=500000, fill_value=None)
absf=np.abs(fcor_nearest)
l_rossby_calc=c_rossby_nearest/absf

btp = nc.variables['BTP'][0,:,:]
btp_nearest = pyresample.kd_tree.resample_nearest(orig_def, btp, \
        targ_def25, radius_of_influence=500000, fill_value=None)

l_rossbyeq_calc=np.sqrt(c_rossby_nearest/(2*btp_nearest))
l_eddy_calc=np.minimum(l_rossby_calc,l_rossbyeq_calc)

### convert regridded radius from cm to km
cm2km=1.e-5
l_eddy_calc_km=l_eddy_calc*cm2km
clevs=[10,20,30,40,50,60,80,100,150,230]
fig=MapZoneContour(targ_def25.lons,targ_def25.lats,l_eddy_calc_km,figsize=(14,9),
                    levels=clevs,
                    llcrnrlon=25.,
                    title="CESM First Baroclinic Rossby Radius (km) $L_{eddy} \sim min(l_r,l_{req})$",
                    fmt='%.1i',
                    imagename='Fig_2_CESM_First_Baroclinic_Rossby_Radius.png')
lrhines = nc.variables['L_RHINES'][0,:,:]
lrhines_nearest = pyresample.kd_tree.resample_nearest(orig_def, lrhines, \
        targ_def25, radius_of_influence=500000, fill_value=None)

###clevs=np.logspace(1,15,40)
###norm=LogNorm()
###fig,ax,cbar =MapContourg(targ_def0,lrhines_nearest,
###                    addzonal=True,
###                    levels1=clevs,
###                    norm=norm,
###                    figsize=(14,5),
###                    title="$L_{rhines} (cmn^2)$ log scale")
###cbar.set_ticks(np.logspace(1,15,15))
if False:
    l_eddy1=np.maximum(l_eddy_calc,lrhines_nearest)
    l_eddy1_km=l_eddy1*cm2km
    clevs=[10,20,30,40,50,60,80,100,150,230]
    fig2=MapZoneContour(targ_def25.lons,targ_def25.lats,l_eddy1_km,figsize=(14,9),
                    levels=clevs,
                    llcrnrlon=25.,
                    title="CESM First Baroclinic Rossby Radius (km) using rhines scale $L_{eddy} \sim max(l_rhines,min(l_r,l_{req}))$",
                    fmt='%.1i' )

In [38]:

Expand Code

In [39]:

Expand Code

$\color{green} {Parameterizing \quad Zonal\quad Eddy\quad Phase\quad Speed}\quad (\color{red}{c})$

NOTES:

$K=u_{rms}∗{\Gamma * L_{eddy} \over (1 + b1 * |u_{mean} - \color{red}{c}|^2 /u_{rms}^2 (z=0)}$
$\require{cancel}\cancel{c = - \beta * {L_r^2}}$ $L_r$ too high at equator
$c = - \beta * L_{eddy}^2$

$\quad\quad\quad\quad\quad$ Eddy Phase speed CESM

In [40]:

Expand Code

# <!-- collapse=True --> 
if False:
    #beta
    clevs=arange(1.e-14,25.e-14,1e-14)
    fig,ax,cbar =MapContourg(targ_def0,btp_nearest,
                    levels1=clevs,
#                    norm=norm,
                    addzonal=True,
                    figsize=(14,5),
                    title="Beta (1/cm*s)",
                    imagename='Fig_3_Beta.png')
l_eddy_calc_sq=l_eddy_calc*l_eddy_calc

if False:
    clevs=np.logspace(1,15,40)
    norm=LogNorm()
    fig,ax,cbar =MapContourg(targ_def0,l_eddy_calc_sq,
                    addzonal=True,
                    levels1=clevs,
                    norm=norm,
                    figsize=(14,5),
                    title="$L_{eddy}^2 (cm^2)$ log scale")
    #                    imagename='Fig_4_L_eddy_SQ.png')

    #cbar.ax.get_xaxis().get_major_formatter().labelOnlyBase = True
    #cbar.ax.get_yaxis().set_major_formatter(plt.LogFormatter(10,  labelOnlyBase=True))
    #ax.yaxis.set_major_locator(ticker.MultipleLocator)
    #ax.yaxis.set_major_formatter(LogFormatter())
    #cbar.formatter.set_scientific(False)
    #cbar.formatter.set_powerlimits((0, 0))
    #cbar.update_ticks()
    cbar.set_ticks(np.logspace(1,15,15))
    fig.savefig('Fig_4_L_eddy_SQ.png',bbox_inches='tight',dpi=200)


c=btp_nearest*l_eddy_calc*l_eddy_calc # units=1/(cm s) * cm^2 = cm/s


clevs=arange(-10,20,1)
fig,ax,cbar = MapContourg(targ_def0,c,
                    addzonal=True,
                    levels1=clevs,
#                    norm=norm,
                    setover='darkred',
                    figsize=(14,5),
                    title="eddy phase speed c = beta*$l_{eddy}^2$ (cm/s) linear scale(NOTE I'm plotting positive c here - changed to negative after this plot")
cbar.set_ticks(np.linspace(-10,20,7))
fig.savefig('Fig_5_eddy_phase_speed.png',bbox_inches='tight',dpi=200)

if False:
    c_rhines=btp_nearest*l_eddy1*l_eddy1 
    clevs=arange(-10,20,1)
    fig,ax,cbar = MapContourg(targ_def0,c_rhines,
                    addzonal=True,
                    levels1=clevs,
#                    norm=norm,
                    setover='darkred',
                    figsize=(14,5),
                    title="eddy phase speed including rhines scale c = beta*$l_{eddy}^2$ (cm/s) linear scale(NOTE I'm plotting positive c here - changed to negative after this plot")
    cbar.set_ticks(np.linspace(-10,20,7))

c=-1.*c

$\quad \quad \quad \quad $Compared to Hughes, the phase speed is much too high near the equator.

$\quad \quad \quad \quad $The first tuning mod is to limit the phase speed to 20cm/s.

$\quad \quad \quad \quad $Eddy Phase speed CESM $(c = max(- \beta * L_{eddy}^2,-20) )$

In [41]:

Expand Code

$\quad\quad\quad\quad\quad\quad$ Hughes Phase Speed (cm/s) from Tulloch Marshall Smith '09

In [42]:

Expand Code

$\color{green}{Zonal\quad Velocity \quad Term}\quad (\color{red}{u_{mean}})$

NOTES:

$K=u_{rms}∗{\Gamma * L_{eddy} \over (1 + b1 * |\color{red}{u_{mean}} - c|^2 /u_{rms}^2 (z=0))}$
Here $U_{mean}$ is just CESM surface velocity.

$\quad\quad \quad \quad\quad $ Zonal Velocity CESM

In [43]:

Expand Code

# <!-- collapse=True --> 
#################  PLOT (u-c)^2  convert m^2/s^2 #######################
umean=nc.variables['U_MEAN'][0,:,:]
umean10_nearest = pyresample.kd_tree.resample_nearest(orig_def, umean,
                targ_def0, radius_of_influence=500000, fill_value=None)

umean10_mps=umean10_nearest*1.e-2

clevs=np.linspace(-60,60,40)
fig,ax,cbar=MapContourg(targ_def0,umean10_nearest,
                    addzonal=True,
                    figsize=(14,5),
                    levels1=clevs,
#                    setover='darkred',
#                    setunder='darkblue',
                    title="U zonal velocity (cm/s)")
#fig.subplots_adjust(left=0.15, right=.9, bottom=0.1, top=0.9);
cbar.set_ticks(arange(-60,60,10))
pylab.xlim([-20.,20.]);
fig.savefig('Fig_7_POP_UVEL.png',bbox_inches='tight',dpi=200)

if False:
    
    # U averaged over 0-50m
    
    clevs=np.linspace(-60,60,40)
    fig,ax,cbar=MapContourg(targ_def0,u50_nearest,
                    addzonal=True,
                    figsize=(14,5),
                    levels1=clevs,
#                    setover='darkred',
#                    setunder='darkblue',
                    title="U 0-50m avg zonal velocity (cm/s)")
    #fig.subplots_adjust(left=0.15, right=.9, bottom=0.1, top=0.9);
    cbar.set_ticks(arange(-60,60,10))
    pylab.xlim([-20.,20.]);




# get 2D versions of the lat and lon variables add longitude start point here!
lon_ecco_orig, lat_ecco_orig = np.meshgrid(lon_ecco[:], lat_ecco[:])

orig_ecco_def = pyresample.geometry.GridDefinition(lons=lon_ecco_orig, lats=lat_ecco_orig)

ecco_u_nearest = pyresample.kd_tree.resample_nearest(orig_ecco_def, ecco_u, \
        targ_def0, radius_of_influence=500000, fill_value=None)
ecco_u_nearest_masked=np.ma.masked_where(ecco_u_nearest<-10, ecco_u_nearest, copy=True)

#print ecco_u_cyc_masked[:,0]

clevs=np.linspace(-60,60,40)
#ecco_u_nearest = pyresample.kd_tree.resample_nearest(orig_def, ecco_u, \
#        targ_def, radius_of_influence=500000, fill_value=None)
fig,ax,cbar=MapContourg(targ_def0,ecco_u_nearest_masked*100,
                    addzonal=True,
                    figsize=(14,5),
                    levels1=clevs,
#                    extend='both',
                    title="ECCO U zonal velocity (cm/s)")
cbar.set_ticks(arange(-60,60,10))
pylab.xlim([-20.,20.]);
fig.savefig('Fig_8_ECCO_UVEL.png',bbox_inches='tight',dpi=200)

$\color{green}{(U-c)\quad Term}$

NOTES:

$c = max(- \beta * L_r^2,-20)$ )

$\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad$ (U-c) CESM

In [44]:

Expand Code

$\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad$ (U-c) Bates et al

In [45]:

Expand Code

$\color{green}{(U-c)^2 \quad Term}$

NOTES:

$c = max(- \beta * L_{eddy}^2,-20)$

$\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad (U-c)^2$ CESM

In [46]:

Expand Code

# <!-- collapse=True --> 

uminusc2sq = nc.variables['U_MINUS_C_SQ'][0,0,:,:]
uminusc2sq_nearest = pyresample.kd_tree.resample_nearest(orig_def, uminusc2sq, \
        targ_def0, radius_of_influence=500000, fill_value=None)
uminusc2sq_mps=uminusc2sq_nearest*1.e-4

if False:
    clevs=np.logspace(-5,0,40)
#    norm=matplotlib.colors.SymLogNorm(linthresh=1e-4, vmin=-1000, vmax=-.001)
    norm=LogNorm()
    fig,ax,cbar=MapContourg(targ_def0,uminusc2sq_mps,
                    addzonal=True,
                    levels1=clevs,
                    norm=norm,
                    figsize=(14,5),
                    title=" ${|(U-c)|}^2$ using new $c_{eddy} {(m^2/s^2)}$ max = %g/min = %g "%(uminusc2sq_mps.max(),uminusc2sq_mps.min()))

    # Customize y tick lables
    yticks=[1e-5,3e-5,1e-4,3e-4,1e-3,3e-3,1e-2,3e-2,1e-1,3e-1,1]
    cbar.set_ticks(yticks)
    cbar.set_ticklabels(yticks)
    pylab.xlim([0.,.1])
    fig.savefig('Fig_12_U_minus_c_SQ_new_log.png',bbox_inches='tight',dpi=200)

clevs=np.linspace(0,.06,25)
fig,ax,cbar=MapContourg(targ_def0,uminusc2sq_mps,
                    addzonal=True,
                    levels1=clevs,
                    setover='darkred',
                    setunder='black',
                    figsize=(14,5),
                    title=" ${|(U-c)|}^2$ using new $c_{eddy} {(m^2/s^2)}$ max = %g/min = %g "%(uminusc2sq_mps.max(),uminusc2sq_mps.min()))
#customize colorbar ticks
yticks=[0,.02,.04,.06]
cbar.set_ticks(yticks)
#customize zonal ticks
pylab.xlim([0,.06])
pylab.xticks([0,.03,.06]);
fig.savefig('Fig_13_U_minus_c_SQ_new.png',bbox_inches='tight',dpi=200)


uminusc2_calc_sq_mps=uminusc2*uminusc2*1.e-4
if False:
    clevs=np.linspace(0,.06,25)
    fig,ax,cbar=MapContourg(targ_def0,uminusc2_calc_sq_mps,
                    addzonal=True,
                    levels1=clevs,
                    setover='darkred',
                    setunder='black',
                    figsize=(14,5),
                    title=" ${|(U-c)|}^2 {(m^2/s^2)}$ calc max = %g/min = %g "%(uminusc2_calc_sq_mps.max(),uminusc2_calc_sq_mps.min()))
    #customize colorbar ticks
    yticks=[0,.02,.04,.06]
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.06])
    pylab.xticks([0,.03,.06]);

$\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad (U-c)^2$ Bates et al

In [47]:

Expand Code

$\color{green}{Parameterizing \quad Eady \quad Growth \quad Rate}\quad \color{red}{\sigma}$

NOTES:

The original derivation of $\sigma_{vi}$ is undefined at the equator because of $f$

$\quad\quad\quad\quad\quad\sigma_{vi} = {f \over \sqrt{R_i}}$

The new derivation of $\sigma_{vi}$:

$\quad\quad\quad\quad\quad R_i = {f^2 N^2 \over {\underbrace{{ {g^2 \over \rho_0^2 } ( \frac{\partial \rho}{\partial y})^2 + {g^2 \over \rho_0^2 } ( \frac{\partial \rho}{\partial z})^2} }_\text{m^4}}}\quad\quad = \quad\quad {f^2N^2 \over m^4} $

$\quad\quad\quad\quad\quad \sigma_{vi} = {f \over \sqrt{R_i}}\quad = \quad {f \over \sqrt{{f^2N^2 \over m^4}}}\quad = \quad {{\cancel f m^2} \over \cancel f N}$

$\quad\quad\quad\quad\quad $Eady Growth Rate $(\sigma) \quad$ CESM

In [48]:

Expand Code

$\quad\quad\quad\quad\quad $Eady Growth Rate $(\sigma) \quad$ Bates et al

In [49]:

Expand Code

$\color{green}{Parameterizing \quad RMS \quad eddy \quad velocity}\quad \color{red}{u_{rms}}$

NOTES:

$K=\color{red}{u_{rms}}∗{\Gamma * L_{eddy} \over (1 + b1 * |u_{mean} - c|^2 /u_{rms}^2(z=0)}$
$u_{rms} = alpha*{\sigma_{vi}}*L_{eddy}$
alpha (scaling constant) = 4
$\sigma_{vi} = {f \over \sqrt{R_i}}\quad = \quad {f \over \sqrt{{f^2N^2 \over m^4}}}\quad = \quad {{\cancel f m^2} \over \cancel f N}$
$u_{rms}$ is limited to 5 cm/s $ max(u_{rms},5.)$

$\quad\quad\quad\quad\quad\quad\quad$ Eddy Velocity $(u_{rms})\quad$ CESM

In [50]:

Expand Code

# <!-- collapse=True --> 
alpha=4
cm2m=.01
#
urms_avg  =  nc.variables['URMS_AVG'][0,:,:]
urms_avg_nearest = pyresample.kd_tree.resample_nearest(orig_def, urms_avg, 
        targ_def0, radius_of_influence=500000, fill_value=None)
urms_avg_mps=urms_avg_nearest*cm2m

if False:
    
    clevs=np.logspace(-3,1,40)
    norm=LogNorm()
    fig,ax,cbar=MapContourg(targ_def0,urms_avg_mps,                
                    addzonal=True,
                    levels1=clevs,
                    norm=norm,
                    figsize=(14,5),
                    setover='darkred',
                    setunder='darkblue',
                    title="$u_{rms} ({m \over s})$ using old sigma log plot(alpha= %.1f,max = %.1f/min = %.2G) "%(alpha,urms_avg_mps.max(),urms_avg_mps.min()))
    #customize colorbar ticks
    yticks=np.logspace(-5,1,7)
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.22])
    pylab.xticks([0,.1,.2]);

    
    
    clevs=np.linspace(0,.4,40)
    fig,ax,cbar=MapContourg(targ_def0,urms_avg_mps,
                    addzonal=True,
                    levels1=clevs,
                    figsize=(14,5),
                    title="$u_{rms} ({m \over s})$ using old sigma  (alpha=4,max = %.1f/min = %.1f) "%(urms_avg_mps.max(),urms_avg_mps.min()))
    #customize colorbar ticks
    yticks=[0,.1,.2,.3,.4]
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.22])
    pylab.xticks([0,.1,.2]);




urms_avg1 = nc.variables['URMS_AVG1'][0,0,:,:]
urms_avg1_nearest = pyresample.kd_tree.resample_nearest(orig_def, urms_avg1, 
        targ_def0, radius_of_influence=500000, fill_value=None)
urms_avg1_mps=urms_avg1_nearest*cm2m

if False:
    clevs=np.logspace(-3,1,40)
    norm=LogNorm()
    fig,ax,cbar=MapContourg(targ_def0,urms_avg1_mps,                
                    addzonal=True,
                    levels1=clevs,
                    norm=norm,
                    figsize=(14,5),
                    setover='darkred',
                    setunder='darkblue',
                    title="$u_{rms} ({m \over s})$ using new sigma log plot(alpha= %.1f,max = %.1f/min = %.2G) "%(alpha,urms_avg1_mps.max(),urms_avg1_mps.min()))
    #customize colorbar ticks
    yticks=np.logspace(-5,1,7)
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.22])
    pylab.xticks([0,.1,.2]);
    fig.savefig('Fig_16_Urms_new_log.png',bbox_inches='tight',dpi=200)

clevs=np.linspace(0,.4,40)
fig,ax,cbar=MapContourg(targ_def0,urms_avg1_mps,
                    addzonal=True,
                    levels1=clevs,
                    figsize=(14,5),
                    title="$u_{rms} ({m \over s})$ using new sigma (alpha= %.1f,max = %.1f/min = %.2G) "%(alpha,urms_avg1_mps.max(),urms_avg1_mps.min()))
#customize colorbar ticks
yticks=[0,.1,.2,.3,.4]
cbar.set_ticks(yticks)
#customize zonal ticks
pylab.xlim([0,.22])
pylab.xticks([0,.1,.2]);
fig.savefig('Fig_17_Urms_new.png',bbox_inches='tight',dpi=200)

cm2m=.01
l_eddy_calc_m=l_eddy_calc*cm2m
urms_calc_mps=alpha*sigma_avg1_nearest*l_eddy_calc_m
if False:
    clevs=np.linspace(0,.4,40)
    fig,ax,cbar=MapContourg(targ_def0,urms_calc_mps,
                    addzonal=True,
                    levels1=clevs,
                    figsize=(14,5),
                    title="$u_{rms} ({m \over s})$ calc value  (alpha= %.1f,max = %.1f/min = %.1f) "%(alpha,urms_calc_mps.max(),urms_calc_mps.min()))
    #customize colorbar ticks
    yticks=[0,.1,.2,.3,.4]
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.22])
    pylab.xticks([0,.1,.2]);

$\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad u_{rms} $ Bates et al

In [51]:

Expand Code

$\color{green}{u_{rms}^2 \quad Term} \quad$

NOTES:

$K={u_{rms}}∗{\Gamma * L_{eddy} \over (1 + b1 * |u_{mean} - c|^2 /\color{red}{u_{rms}^2} (z=0)} $
$u_{rms}^2$ is a surface value

$\quad\quad\quad\quad\quad\quad\quad$ Eddy Velocity Squared $(u_{rms}^2)\quad$ CESM

In [52]:

Expand Code

#### <!-- collapse=True --> 
urms_avg1_sq = nc.variables['URMS_AVG1_SQ'][0,0,:,:]
urms_avg1_sq_nearest = pyresample.kd_tree.resample_nearest(orig_def, urms_avg1_sq, 
        targ_def0, radius_of_influence=500000, fill_value=None)
urms_avg1_sq_mps=urms_avg1_sq_nearest*1.e-4
###urms_avg1_sq_mps =urms_avg1_mps*urms_avg1_mps

urms_avg_sq = nc.variables['URMS_AVG_SQ'][0,:,:]
urms_avg_sq_nearest = pyresample.kd_tree.resample_nearest(orig_def, urms_avg_sq, 
        targ_def0, radius_of_influence=500000, fill_value=None)
urms_avg_sq_mps=urms_avg_sq_nearest*1.e-4

if False:
    clevs=np.logspace(-5,1,42)
    norm=LogNorm()
    fig,ax,cbar=MapContourg(targ_def0,urms_avg_sq_mps,                
                    addzonal=True,
                    levels1=clevs,
                    norm=norm,
                    figsize=(14,5),
                    setover='darkred',
                    setunder='darkblue',
                    title="$u_{rms}^2 old sigma({m^2 \over s^2})$ log plot (alpha=%.1f,max = %.2G/min = %.2G) "%(alpha,urms_avg_sq_mps.max(),urms_avg_sq_mps.min()))
    #customize colorbar ticks
    yticks=np.logspace(-5,1,7)
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.06])
    pylab.xticks([0,.03,.06]);



    ##clevs=np.linspace(0,.06,40)
    clevs=np.linspace(.001,.06,40)
    fig,ax,cbar=MapContourg(targ_def0,urms_avg_sq_mps,
                    addzonal=True,
                    levels1=clevs,
                    setover='darkred',
                    setunder='darkblue',
                    figsize=(14,5),
                    title="$u_{rms}^2 old sigma ({m^2 \over s^2})$ (alpha=%.1f,max = %.2G/min = %.2G) "%(alpha,urms_avg_sq_mps.max(),urms_avg_sq_mps.min()))
    #customize colorbar ticks
    yticks=[.001,.02,.04,.06]
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.06])
    pylab.xticks([0,.03,.06]);


    
    
    clevs=np.logspace(-5,1,42)
    norm=LogNorm()
    fig,ax,cbar=MapContourg(targ_def0,urms_avg1_sq_mps,                
                    addzonal=True,
                    levels1=clevs,
                    norm=norm,
                    figsize=(14,5),
                    setover='darkred',
                    setunder='darkblue',
                    title="$u_{rms}^2 new sigma({m^2 \over s^2})$ log plot (alpha=%.1f,max = %.2G/min = %.2G) "%(alpha,urms_avg1_sq_mps.max(),urms_avg1_sq_mps.min()))
    #customize colorbar ticks
    yticks=np.logspace(-5,1,7)
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.06])
    pylab.xticks([0,.03,.06]);
    fig.savefig('Fig_18_Urms_SQ_log.png',bbox_inches='tight',dpi=200)



##clevs=np.linspace(0,.06,40)
clevs=np.linspace(.001,.06,40)
fig,ax,cbar=MapContourg(targ_def0,urms_avg1_sq_mps,
                    addzonal=True,
                    levels1=clevs,
                    setover='darkred',
                    setunder='darkblue',
                    figsize=(14,5),
                    title="$u_{rms}^2 new sigma ({m^2 \over s^2})$ (alpha=%.1f,max = %.2G/min = %.2G) "%(alpha,urms_avg1_sq_mps.max(),urms_avg1_sq_mps.min()))
#customize colorbar ticks
yticks=[.001,.02,.04,.06]
cbar.set_ticks(yticks)
#customize zonal ticks
pylab.xlim([0,.06])
pylab.xticks([0,.03,.06]);
fig.savefig('Fig_19_Urms_SQ.png',bbox_inches='tight',dpi=200)

urms_calc_mps_sq=urms_calc_mps*urms_calc_mps

if False:
    clevs=np.linspace(.001,.06,40)
    fig,ax,cbar=MapContourg(targ_def0,urms_calc_mps_sq,
                    addzonal=True,
                    levels1=clevs,
                    setover='darkred',
                    setunder='darkblue',
                    figsize=(14,5),
                    title="$u_{rms}^2 ({m^2 \over s^2})$ (alpha=%.1f,max = %.2G/min = %.2G) "%(alpha,urms_calc_mps_sq.max(),urms_calc_mps_sq.min()))
    #customize colorbar ticks
    yticks=[.001,.02,.04,.06]
    cbar.set_ticks(yticks)
    #customize zonal ticks
    pylab.xlim([0,.06])
    pylab.xticks([0,.03,.06]);

$\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad u_{rms}^2 $ Bates et al

In [53]:

Expand Code

$\color{green}{{(U-c)}^2 \over {u_{rms}^2}} \quad \color{green}{Term}$

NOTES:

$K={u_{rms}}∗{\Gamma * L_{eddy} \over (1 + b1 * \color{red}{|u_{mean} - c|^2/{u_{rms}^2} (z=0)}} $
$c = max(- \beta * L_{eddy}^2,-20)$
$u_{rms} = max(u_{rms},5.)$
$u_{rms}^2$ is a surface value

$\quad\quad\quad\quad\quad\quad$ ${(U-c)}^2 \over {u_{rms}^2} \quad$ CESM

In [54]:

Expand Code

$\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad {|u-c|^2 \over u_{rms}^2} $ Bates et al

In [55]:

Expand Code

$\color{green}{\quad Suppression \quad factor = {1 \over (1 + b1 * |\bar u - c|^2 /u_{rms (z=0)}^2 )}}$

NOTES:

$c = max(- \beta * L_{eddy}^2,-20)$
$u_{rms} = max(u_{rms},5.)$
$u_{rms}^2$ is a surface value
$b1$ (scaling constant) = 4.

$\quad\quad\quad\quad\quad\quad\quad$ Suppression Factor CESM

In [56]:

Expand Code

$\quad\quad\quad\quad\quad\quad\quad$ Suppression Factor Bates

In [57]:

Expand Code

$\quad\quad\quad\quad\quad\quad\quad$ Suppression Factor (Zonal x depth) CESM

In [58]:

Expand Code

clevs=np.linspace(0,1,11)

z_t= ncregrid.variables['z_t'][:]

latregrid=ncregrid.variables['lat'][:]
lonregrid=ncregrid.variables['lon'][:]
lon=lonregrid

suppregrid= ncregrid.variables['SUPP'][0,:,:,:]
suppregridz=suppregrid.mean(axis=2)
clevs=np.linspace(0,1,100,endpoint=True)
yticks=[0,.2,.4,.6,.8,1.]
fig,ax,cbar=plotzonalocn(latregrid, z_t, suppregridz,
                    levels1=clevs,
                    levels2=clevs,
#                    equalcontourscale=True,
                    figsize=(10.5,6),
                    title="Suppression factor (custom contours max = %.1f/min = %.1f) "%(suppregridz.max(),suppregridz.min()))

# Add axis labels and plot title                                                                                                                                                                                              
plt.xlabel('Latitude')                                                                                                                                                                                                
plt.ylabel('Depth (cm)')                                                                                                                                                                                              

# Customize y tick lables                                                                                                                                                                                                                                                                                                                                                                                                              
cbar.set_ticks(yticks)                                                                                                                                                                                                
cbar.set_ticklabels(yticks)                                                                                                                                                                                           
cbar.set_label('fraction')                                                                                                                                                                                            
fig.savefig('Fig_23_Suppression_zonal_depth.png',bbox_inches='tight',dpi=200)

$\quad\quad\quad\quad\quad\quad\quad$ Suppression Factor (Zonal x depth) Bates

In [59]:

Expand Code

$\color{green} {\quad Parameterizing \quad The \quad Mixing \quad Term \quad }\color{red}{L_{mix}}$

NOTES:

$\color{red}{L_{mix}} = \Gamma * L_{eddy} * Suppression$
$Suppression= {1 \over (1 + b1 * |u_{mean} - c|^2 /u_{rms (z=0)}^2 )}$
$\color{red}{L_{mix}} = {\Gamma * L_{eddy} \over (1 + b1 * |u_{mean} - c|^2 /u_{rms}^2 (z=0)}$
$\Gamma = 1.75$ (Tuning mod: the original Gamma of .35 produced a Kappa with the correct structure but too weak.

$\quad\quad\quad\quad\quad\quad\quad$ LMIX CESM

In [60]:

Expand Code

$\color{green}{\quad Parameterizing \quad Eddy \quad Diffusivity \quad} (\color{red}{K})$

NOTES:

$\color{red}{K}=u_{rms}∗L_{mix}$

$\quad\quad\quad\quad\quad\quad\quad$ Eddy Diffusivity (K) CESM

In [61]:

Expand Code

$\quad\quad\quad\quad\quad\quad\quad$ Eddy Diffusivity (K) Bates

In [62]:

Expand Code

$\quad\quad\quad\quad\quad\quad\quad$ (Zonal x Depth) CESM

$\quad\quad\quad\quad\quad\quad\quad\quad\quad$ K limited to (100 < K < 10000)

In [63]:

Expand Code

# <!-- collapse=True --> 

###ncregrid = Dataset(basedir+'/steering/g.e11.GIAF.T62_gx1v6.steer.017.regrid.nc')
z_t= ncregrid.variables['z_t'][:]

latregrid=ncregrid.variables['lat'][:]
lonregrid=ncregrid.variables['lon'][:]
lon=lonregrid

kvert= ncregrid.variables['KVERT'][0,:,:,:]
kvertz=kvert.mean(axis=2)/1.e4

clevs=np.arange(0,100,20)
clevs=np.append(clevs,np.arange(100,1000,100))
clevs=np.append(clevs,np.arange(1000,4000,1000))
fig,ax,cbar=plotzonalocn(latregrid, z_t, kvertz,
                    levels1=clevs,
                    levels2=clevs,
                    equalcontourscale=True,
                    figsize=(10.5,6),
                    title="Eddy diffusivity $K$ (custom contours max = %.1f/min = %.1f) "%(kvertz.max(),kvertz.min()))

# Add axis labels and plot title                                                                                                                                                                                              
plt.xlabel('Latitude')                                                                                                                                                                                                
plt.ylabel('Depth (cm)')                                                                                                                                                                                              

# Customize y tick lables                                                                                                                                                                                                     
yticks=clevs                                                                                                                                                                                                          
cbar.set_ticks(yticks)                                                                                                                                                                                                
cbar.set_ticklabels(yticks)                                                                                                                                                                                           
cbar.set_label('m**2/s')                                                                                                                                                                                            
fig.savefig('Fig_29_K_zonal_depth_custom.png',bbox_inches='tight',dpi=200)

clevs=np.linspace(0,2500,100,endpoint=True)
yticks=[0,500,1000,1500,2000,2500]
fig,ax,cbar=plotzonalocn(latregrid, z_t, kvertz,
                    levels1=clevs,
                    levels2=clevs,
#                    equalcontourscale=True,
                    figsize=(10.5,6),
                    title="Eddy diffusivity $K$ (Linear Scale max = %.1f/min = %.1f) "%(kvertz.max(),kvertz.min()))

# Add axis labels and plot title                            
levs=[500,1000,1500]
cs=ax.contour(latregrid, z_t, kvertz,levels=levs,colors='w',linewidths=2)
plt.xlabel('Latitude')                                                                                                                                                                                                
plt.ylabel('Depth (cm)')                                                                                                                                                                                              
plt.clabel(cs, inline=1, fontsize=15)
# Customize y tick lables                                                                                                                                                                                                                                                                                                                                                                                                            
cbar.set_ticks(yticks)                                                                                                                                                                                                
cbar.set_ticklabels(yticks)                                                                                                                                                                                           
cbar.set_label('m**2/s')                                                                                                                                                                                            
fig.savefig('Fig_30_K_zonal_depth_linear.png',bbox_inches='tight',dpi=200)

Here is the Bate's version

In [64]:

Expand Code

$\quad\quad\quad\quad\quad\quad\quad$ Zonal average of the N2 normalized scaling CESM

In [65]:

Expand Code

clevs=np.linspace(0,1,11)
###ncregrid = Dataset(basedir+'/steering/g.e11.GIAF.T62_gx1v6.steer.015.regrid.nc')
n2norm=ncregrid.variables['BFREQ_SQ_NORM'][0,:,:,:]
n2normz=n2norm.mean(axis=2)
fig,ax,cbar=plotzonalocn(latregrid, z_t, n2normz,
                    levels1=clevs,
                    levels2=clevs,
                    equalcontourscale=True,
                    figsize=(14,5),
                    title="N2 Scaling (max = %.1f/min = %.1f) "%(n2normz.max(),n2normz.min()))

# Add axis labels and plot title                                                                                                                                                                                              
plt.xlabel('Latitude')                                                                                                                                                                                                
plt.ylabel('Scaling factor')                                                                                                                                                                                              

# Customize y tick lables                                                                                                                                                                                                     
yticks=clevs                                                                                                                                                                                                          
cbar.set_ticks(yticks)                                                                                                                                                                                                
cbar.set_ticklabels(yticks)                                                                                                                                                                                           
cbar.set_label('fraction')                                                                                                                                                                                            

fig.savefig('Fig_31_N2_NORM_zonal_depth.png',bbox_inches='tight',dpi=200)

In [66]:

Expand Code

Out[66]:

The alternate derivation of $\sigma_{vi}$:

$ R_i = {N^2 }$

$ N^2 = {{-g _0 }} $

After hydrostatic and geostrophic approximations

$f \frac {\partial v}{\partial z} = {{-g \over \rho_0 }\frac {\partial \rho}{\partial x}}; \quad\quad f \frac {\partial u}{\partial z} = {{g \over \rho_0 }\frac {\partial \rho}{\partial y}} $

so

$\frac {\partial v}{\partial z} = {{-1\over f}{g \over \rho_0 }\frac {\partial \rho}{\partial x}}; \quad\quad \frac {\partial u}{\partial z} = {{1\over f}{g \over \rho_0 }\frac {\partial \rho}{\partial y}}$

$R_i = {f^2 N^2 \over {\underbrace{{ {g^2 \over \rho_0^2 } ( \frac{\partial \rho}{\partial y})^2 + {g^2 \over \rho_0^2 } ( \frac{\partial \rho}{\partial z})^2} }_\text{m^4}}}\quad\quad = \quad\quad {f^2N^2 \over m^4} $

$\sigma_{vi} = {f \over \sqrt{R_i}}\quad = \quad {f \over \sqrt{{f^2N^2 \over m^4}}}\quad = \quad {{\cancel f m^2} \over \cancel f N}$

$RX_1 = RX_{east} = \Delta\rho_x = \rho_{i+1,j} - \rho_{i,j}$

$RY_1 = RY_{north} = \Delta\rho_y = \rho_{i,j+1} - \rho_{i,j}$

$RZ_1 = RZ_{k+1} = \Delta\rho_z = \rho_{k} - \rho_{k+1}$

$\displaystyle{1 \over L_{R_i}} \displaystyle\int_{2000m}^{100m} \left\lbrace { {-g\over\rho_0}{\frac {\partial \rho} {\partial z}} \over { {g^2 \over \rho_0^2 } \left[( \frac{\partial \rho}{\partial y})^2 + ( \frac{\partial \rho}{\partial z})^2\right] } } \right\rbrace dz$

Note: missing $f^2$ which will be cancelled when forming $\sigma_{vi}$

$\quad\quad$ so $\cdots$ this is not $R_i$

Implementation notes

Numerator : Top $= -grav * RZ_{SAVE}(\cdots k+1) * dzwr(k)$

Denominator :

$\begin{align} work1 = p25 & * ( RX(..,i_{east},k)^2 \\ & + RX(..,i_{west},k)^2 \\ & + RX(..,i_{east},k+1)^2 \\ & + RX(..,i_{west},k+1)^2 ) / DXT(i,j)^2 \\ \end{align}$

$\begin{align} work2 = p25 & * ( RY(..,j_{north},k)^2 \\ & + RY(..,j_{south},k)^2 \\ & + RY(..,j_{north},k+1)^2 \\ & + RY(..,j_{south},k+1)^2 / DYT(i,j)^2 \\ \end{align}$

$\begin{align} work3 = {\left( TOP \over (grav^2*(work1+work2))\right)}*dzw(k) \end{align}$

Notes:

1)Need to be careful at top and bottom of ocean
2)Accurate dzw(k) for each (i,j) to form $L_{R_i}$
3) When constructing $sigma$ itself, use $RZ_{SAVE}$ with a minimum N value
4) use eps2

In []:

POP Implementation of Mesoscale eddy diffusivity based on mixing theory

Development Status

POP standard diagnostics for the new steering case

Implementation Notes on Diffusivity with Steering Level Suppression (Gokhan and John Marshall)

\(\quad\quad\quad\quad\quad\quad\) Eddy Length Scales:

\(\quad\quad\quad\quad\quad\quad\) Eddy Velocity:

\(\quad\quad\quad\quad\quad\quad\) Zonal phase speed:

References:

Questions:

\(\color{green}{\quad CESM \quad Implementation \quad of \quad Steering \quad Parameterization}\)

\(\color{green} {Parameterizing \quad the\quad Eddy\quad Length\quad Scale} \)

\(\quad\quad\quad\quad\quad\quad\quad\) First Baroclinic Wave Speed \(c_r\) CESM

\(\quad\quad\quad\quad\quad\quad\) Chelton Baroclinic gravity wave phase speed

\(\quad\quad\quad\quad\quad\quad\) First Baroclinic Rossby Radius CESM \(\quad (L_{eddy})\)

\(\color{green} {Parameterizing \quad Zonal\quad Eddy\quad Phase\quad Speed}\quad (\color{red}{c})\)

\(\quad\quad\quad\quad\quad\) Eddy Phase speed CESM

\(\quad \quad \quad \quad \)Eddy Phase speed CESM \((c = max(- \beta * L_{eddy}^2,-20) )\)

\(\quad\quad\quad\quad\quad\quad\) Hughes Phase Speed (cm/s) from Tulloch Marshall Smith '09

\(\color{green}{Zonal\quad Velocity \quad Term}\quad (\color{red}{u_{mean}})\)

\(\quad\quad \quad \quad\quad \) Zonal Velocity CESM

\(\color{green}{(U-c)\quad Term}\)

\(\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\) (U-c) CESM

\(\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\) (U-c) Bates et al

\(\color{green}{(U-c)^2 \quad Term}\)

\(\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad (U-c)^2\) CESM

\(\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad (U-c)^2\) Bates et al

\(\color{green}{Parameterizing \quad Eady \quad Growth \quad Rate}\quad \color{red}{\sigma}\)

\(\quad\quad\quad\quad\quad \)Eady Growth Rate \((\sigma) \quad\) CESM

\(\quad\quad\quad\quad\quad \)Eady Growth Rate \((\sigma) \quad\) Bates et al

\(\color{green}{Parameterizing \quad RMS \quad eddy \quad velocity}\quad \color{red}{u_{rms}}\)

\(\quad\quad\quad\quad\quad\quad\quad\) Eddy Velocity \((u_{rms})\quad\) CESM

\(\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad u_{rms} \) Bates et al

\(\color{green}{u_{rms}^2 \quad Term} \quad\)

\(\quad\quad\quad\quad\quad\quad\quad\) Eddy Velocity Squared \((u_{rms}^2)\quad\) CESM

\(\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad u_{rms}^2 \) Bates et al

\(\color{green}{{(U-c)}^2 \over {u_{rms}^2}} \quad \color{green}{Term}\)

\(\quad\quad\quad\quad\quad\quad\) \({(U-c)}^2 \over {u_{rms}^2} \quad\) CESM

\(\quad\quad\quad\quad\quad\quad\quad\quad\quad\quad {|u-c|^2 \over u_{rms}^2} \) Bates et al

\(\color{green}{\quad Suppression \quad factor = {1 \over (1 + b1 * |\bar u - c|^2 /u_{rms (z=0)}^2 )}}\)

\(\quad\quad\quad\quad\quad\quad\quad\) Suppression Factor CESM

\(\quad\quad\quad\quad\quad\quad\quad\) Suppression Factor Bates

\(\quad\quad\quad\quad\quad\quad\quad\) Suppression Factor (Zonal x depth) CESM

\(\quad\quad\quad\quad\quad\quad\quad\) Suppression Factor (Zonal x depth) Bates

\(\color{green} {\quad Parameterizing \quad The \quad Mixing \quad Term \quad }\color{red}{L_{mix}}\)

\(\quad\quad\quad\quad\quad\quad\quad\) LMIX CESM

\(\color{green}{\quad Parameterizing \quad Eddy \quad Diffusivity \quad} (\color{red}{K})\)

\(\quad\quad\quad\quad\quad\quad\quad\) Eddy Diffusivity (K) CESM

\(\quad\quad\quad\quad\quad\quad\quad\) Eddy Diffusivity (K) Bates

\(\quad\quad\quad\quad\quad\quad\quad\) (Zonal x Depth) CESM

Here is the Bate's version

\(\quad\quad\quad\quad\quad\quad\quad\) Zonal average of the N2 normalized scaling CESM

The alternate derivation of \(\sigma_{vi}\):

Implementation notes